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We study two-dimensional quantum dots using the variational quantum Monte Cario technique 
in the weak-confinement limit where the system approaches the Wigner molecule, i.e., the classical 
solution of point charges in an external potential. We observe the spin-polarization of electrons 
foUowed by a smooth transition to a Wigner-molecule-like state as the confining potential is made 
weaker. 
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I. INTRODUCTION 

Semiconductor quantum dots (QD) are small devices 
containing a tunable njjmber of electrons in an external 
confinement potential.EI The significant progress in the. 
fabrication of these devices during the last few yearstl 
has stimulated an increasing interest in investigating the 
properties of such systems. From the theoretical point 
of view, QDs are ideal many-electron objects for the 
study of fundamental physical properties of correlated 
electrons. 

Perhaps the most striking feature in these artificial 
atoms is that the system parameters can easily be 
changed, unlike in real atoms where the parameters are 
natural constants. Also the typical length scales of inter- 
actions and confinement are equal, which should in prin- 
cipie make the correlation effects more enhanced com- 
pared to normal atoms. However, the experimentally ob- 
served states of QDs in weak magnetic fields can easily 
be understood as sinsfle configurations of non-interacting 
one-particle statesBtl One might suppose that the lack of 
a strong central Coulomb potential (which organizes the 
states of normal atoms to shells) would make it more diffi- 
cult to identify the electronic structure of artificial atoms. 
It turns out, however, that the important role of the cen- 
tral potential is takejx. by the symmetries of the many- 
body wave function.Ela Consequently the understanding 
of the topology of the many-body wave function is of 
central importance. 

A most suitable method for studying the many-body 
wave functions is the variational quantum Monte Cario 
(VMC) technique. It is intimately coupled to the struc- 
ture of the wave function, and the use of the variational 
principie enables one to find the best form. The exact 
diagonalization (ED) method is useful for the studies of 
small electrón numbersci or for finding eigenstates in some 
important subspace of the fuU Hilbert space, as in the 
studies of electrons in such a strong magnetic field that 
the physics is mainly determined. by the single-particle 
states of the lowest Landau lévela ED has been used for 
QDs by many authors, see Ref. || and references therein. 
Mean-field methods such as the spin-density functional 
theory (SDFT) axje espccially useful when the number of 
particles is large.El 

In this work, we first compare our VMC results with 



the most accurate energies in the literature, showing that 
a rather intuitive wave function results in extremely accu- 
rate energies. Then we concéntrate on the weak confine- 
ment limit of a six-electron QD. One should note that 
making confinement weaker makes the density smaller. 
In this limit the QD approaches the classical solution 
of point charges in the QD potential. The most accu- 
rate results |£or this limit are restricted to rather strong 
confinemeniü, where no trace of Wigner molecule or 
spin-polarization was found. For weaker confinements 
than considered in Ref. |5| we find that the system spin- 
polarizes before the electrons localize around the classi- 
cal positions. We also introduce a conditional probability 
density which shows to be a good measure of the extent 
of the localization of the electrons. In addition, we find 
an interesting independent-electron-type scaling of the 
energy for a large range of confinement strengths. 



II. MODEL AND METHOD 

The commonly used Hamiltonian for N electrons in a 
two-dimensional QD can be written as 
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where w is the strength of the external confining poten- 
tial, and the effective mass of the electrons m* and the 
dielectric constant e are used to model the properties of 
the semiconductor material studied. For GaAs, the,ma- 
terial parameters are m* /m^ — 0.067 and e = 12.4J3 We 
assume that the electrons move in the z = O plañe and 
omit the magnetic field in this study. 

In this work, we use variational wave functions of the 
form 
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where the two first factors are Slater determinants for 
the two spin types, and J is a Jastrow two-body cor- 
relation factor. We neglect the three-body and higher 
correlations. For the Jastrow factor we use 
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where a is fixed by the cusp condition to be 3 a for pair 
of equal spins and 1 for opposite ones and 6 is a parame- 
ter, difFerent for both spin-pair possibilities. We take the 
wave function to be real, as we have neglected the mag- 
netic field. The single-particle states tp are expanded in 
the basis of gaussians as 

i;{r) = J2 (^)^««. (y)^'^'' ' (4) 

i=l 

where H„ is a Hermite polynomial of order n and f — 
(i, y) = {cxXjCyn), where Cx and Cy scale coordinates. 

We calcúlate the energy E which is bound by the exact 
energy Eq using 

1 *^ 

Eo<E= tjE^l(RO , (5) 
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where the local energy is í?l = H'í'/'íi, and depends on 
the electrón configuration R. The configurations R are 
distributed according to jí*!^. We optimize the varia- 
tional parameters using the stochastic gradient approxi- 
mation (SGA).EI The SGA optimization method involves 
stochastic simulation in two spaces: the configuration 
and the parameter space. In the configuration space, a 
set of m configurations {Rj} is sampled from a distribu- 
tion |Í'(q:)P, where a is the current parameter vector. 
In the parameter space, the parameters at iteration i + 1 
are obtained from the previous ones by the formula 

Oíi+i ^ oti~ %Vo¿Qi , (6) 

where % is a scaling factor of the step length and Q 
is an approximation to the cost function. For energy 
minimization the cost function is simply the mean of the 
local energies over the set of configurations: 

2 = (^l) = -E^l(R,) . (7) 

The scaling factor 7 has an important role in averaging 
out the fluctuations in the approximate gradient, ensur- 
ing the convergence. On the othcr hand, too small a 
valué of 7 would overdamp the simulation. If one uses a 
sequence 'jiJX , one should have ^ < /3 < 1. 

Lin et alx3 have shown that in the case of real wave 
functions and energy minimization, the derivative of the 
energy E with respect to a variational parameter a¿ is 
simply 

a^r 'H^ [H^ n • 

where the ^average (...) is over the whole Metrópolis 
simulation.Ej One can implement this simple formula also 
for the SGA algorithm, with the small modification that 
the average is taken over only the current set of m con- 
figurations. 



III. RESULTS 

A. Comparison with other approaches 

Before presenting the results for the Wigner-molecule- 
limit of a six-electron QD, we compare our VMC results 
with the most accurate QD results up-to-date. In a re- 
cent yMC and diffusion quantum Monte Garlo (DMC) 
study,Ll Pederiva et al. study quantum dots using a sim- 
ilar model as we do. They find the accuracy of VMC to 
be rather good compared to DMC, except in the case 
of three electrons. For this case, we have first taken 
the single-particle states to be the non-interacting ones 
with quantum numbers (0,0) and (1,0). With the GaAs 
cholee of System parameters given above, the confine- 
mept energy corresponding to the study of Pederiva et 
al.U is around 3_32 meV. The total energy is found to be 
26.563(1) meWJ which is reasonably cióse to the DMG 
valué of Pederiva et ai, namely 26.488(3) meV.Q On the 
other hand, the VMG energy reported by Pederiva et al. 
is 29.669(3) meV.IZl Optimizing the exponentials lowers 
our VMG energy to 26.5406(8) meV. The difference of 
our VMC energy to the DMG one is around 0.05 meV 
which is small compared to the SDFT error ~ 0.4 meV.Q 
For the six-electron case, the energy with non-interacting 
single-particle states is found to be 90.27(1) meV, which 
is again clo,ser to the DMC energy 90.11(1) meV of Ped- 
eriva et aM than the VMC one 90.368(4) meV. By free- 
ing the parameters in the single particle states one does 
not lower the energy within the statistical error, and the 
optimal valúes are thus equal to the ones with the non- 
interacting states. This is a very important result, show- 
ing that the change in the wave function introduced by 
the Coulomb interaction is very accurately taken into ac- 
count by the two-body Jastrow factors used. The reason 
why the optimization of the single-particle states was im- 
portant in the three-electron case is most probably that 
there the number of spin-up and spin-down electrons is 
different. 

It is also interesting to compare the results obtained 
with VMG ■mth those of Reimann et al. for the six- 
electron case.El In their study, they use both SDFT and 
ED and consider six electrons with various strengths 
of the confining potential. We compare four different 
strengths of confinement, corresponding in their work to 
the cases of =1, 2, 3, and 4 a¿. Their ED energies 
are given in Table || with our VMC energies. One can 
see that the difference in the energies is between 0.1 and 
0.3 meV, the VMC energies being lower. This compari- 
son shows that the finite basis used by Reimann et al. is 
too restricted to describe the many-body wave function 
accurately. This comparison also shows that the results 
obtained with VMC are very accurate for the particle 
number in question. 

We also compare the accuracy of our VMC results with 
the path-intesral Monte Cario simulations of Reusch et 
al. for iV = 8.I3 As we are below mainly interested in the 
fuUy spin-polarized states, we compare only the S = A 
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TABLE I: Total energy (in meV) of the six-electron dot for 
difíerent confinements and the two spin states. 





S=0 




S=3 






VMC 


ED" 


VMC 


ED" 


7.576 


168.90(1) 


169.2 


180.40(1) 


180.5 


2.678 


76.91(1) 


77.17 


79.271(4) 


79.38 


1.458 


49.101(5) 


49.35 


49.934(3) 


50.10 


0.947 


35.864(3) 


36.11 


36.231(2) 


36.41 
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energies. This is not a closed-shell case, and our vari- 
ational wave function for the two highest states is con- 
structed from the four states with Ux + Uy = 3. For 
interaction strength C = 2, their energy is 48.3(2) Hlo 
which is shghtly higher (but within error bars) than our 
energy of 48.201(1) huo. As one can see, the statistical 
error is two orders of magnitude smaller in our result. 
For their most strongly interacting case, namely C = 8, 
Reusch eí al. obtain energy 103.26(5) hu¡ which is again 
less accurate than our energy 103.137(1) fiu¡. Also this 
test shows that our results are very accurate in the hmit 
of strong interaction. 

One can conclude from these comparisons that the 
wave function used is very accurate. Most striking is 
the observation that after the Jastrow factor is added to 
the wave function, the remaining problem can be treated 
using non-interacting single-particle states. Usually in 
VMC studies this remaining problem contains some in- 
teractions which are taken into account in a mean-field 
fashion by modified single-particle states. 



B. Classical limit 

Next, we study the transition to the classical limit 
in the six-electron case. The ground-state structure of 
six purely classical point charges in a parabolic potential 
minimizes the energy 
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where we have used reduced unitsH. If we keep other pa- 
rameters fixed and change only the confinement strength 
Lü (as we do below), one can see that the interaction 
strength C scales as C oc o;"^/^. The minimum-energy 
positions of electrons form a pentagon around one elec- 
trón at the center. One can find the scaling of the classi- 
cal cluster size and energy by writing the coordínales as 
r = Tcf , where coordinates f are fixed and the scaling is 
in Te- This results for the energy 
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where Vi and V2 are constants, and solving dEch/drc = 
O results in = {CV2/Vif/^ oc w-^/^. Thus the en- 
ergy scales (in units of lo) as Ecl oc uj" 
imum energy E^^ can be found to be 
30.46(?iw[meV])2/3. 

The quantum-mechanical energies resulting from our 
VMC calculations for two spin polarizations are pre- 
sented in Fig. |l| with the classical energy E^^. One can 
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FIG. 1: Total energy for spin states S = ?> (marked with 
pluses) and 5 = (circles) as a function of Tiu>. The line 
presents the classical energy -Eol- 

see that the quantum-mechanical energies are very cióse 
to each other, especially in the small Tioj limit, where 
the energies also approach the classical one. The differ- 
ence between the quantum-mechanical energies can be 
seen more clearly in Fig. H, where we have plotted the 




0.4 0.6 
ñco [meV] 

FIG. 2: Energy difFerence between the spin states S = 3 and 
5 = O as a function of hu. The line is to guide the eye. 

energy difference between the fuUy- and non-polarized 
states. One can see that spin-polarization is predicted at 
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Tiuj w 0.28 meV. We have not found ground states with 
partial spin-polarization. Below, we concéntrate on the 
fully spin-polarized case. If one subtracts from the total 
energy the minimum valué of the classical potential en- 
ergy E^j^ , the rcmaining energy has a linear behavior at 
small Tiüj as shown in Fig. ^. One way of understanding 



probability distribution p(r), defined as 
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FIG. 3: Remaining energy for S = 3 when -EJi (see text) 
is subtracted. The Une presents a linear fit to points ñw < 
0.1 meV. 



this linear behavior is to consider first the limit hto 0. 
In this limit, the electrons localizo to the classical posi- 
tions. If hüj is now made largor, the electrons start to 
oscillate around the classical positions. A first approx- 
imation for this is to assume that each electrón is in a 
harmonio potential V{ri) = ^Cj'^rf, with the strength 
scaling as w oc w as a function of uj. As a result of this, 
the total energy has, apart from the classical energy, a 
zero-point energy that also scales as w. Thus we obtain 
a linear behavior for the difference between the total en- 
ergy and the classical energy similarly as in Fig. ^. The 
error in the linear fit of Fig. || is smaller than 3 ficY for 
huj < 0.1 meV where the fit is made, and for huj = 1 nieV 
still only 0.2 meV, which is of the same order as the dif- 
ference in energy between the VMC and ED results in 
Table||. For hw = 2 meV the error is around 1 meV. Thus 
the approximation derived above works surprisingly well 
even for a rather strong confinement. One should note 
that the approximation did not contain any Information 
of the interaction between electrons (apart from the clas- 
sical potential energy), and thus the system can be seen 
energetically as a collection of nearly independent elec- 
trons oscillating in an effective potential. 

The most probable configuration R*, maximizing the 
density |Í'(R)|^ should approach in the limit of weak 
confinement the classical electrón positions. This is not, 
however, enough to show that the system is cióse to a 
classical one. One can study the quantum fluctuations 
very conveniently using the conditional single-particle 
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where the coordinates r* are fixed to the ones from the 
most probable configuration R*. In the classical limit, 
the density p(r) is more and more peaked around the clas- 
sical position r*, but still shows quantum fluctuations. 
For the two-electron case, this is very similar to the con- 
ditiouíil probability distribution used for a two-electron 
QD.E3 One difference is that in p the fixed electrón is 
at the most probable position, which is in our opinión 
the most natural cholee. The most important advan- 
tage of p is that it shows much more clearly the amount 
of localization for larger partióle numbers than the con- 
ditional probability distribution. The reason for this is 
that one usually fixes only one electrón in constructing 
the conditional probability distribution, and when the 
partióle number is large, the effect of one fixed electrón 
gets smaller, and the rest of electrons can, for example, 
show collective motion that conceals the localization. If, 
on the other hand, one fixes all but one electrón, the 
most natural cholee is p. The calculation of p is very 
easy, especially in VMC. One should first, of course, find 
the most probable electrón positions. In doing this, the 
gradient of the wave function (needed also for the cal- 
culation of the local energy, and for this reason usually 
done analytically) is very useful. After that, one moves 
the "probé electrón" to all points where the valué of p is 
wanted, and evalúales the rallo of wave funclions as in 
Eq. ( |Tl| ) . This rallo is automatically done when sampling 
the configurations in a VMC simulation. One should also 
notice that p does not contain noise unlike many more 
common VMC observables, such as the density or the 
radial pair distribution function. 

In Fig. ^ we show p for four different confinement 
strenglhs from 10 to 0.01 mcV. The most probable elec- 
trón positions are for all strenglhs similar to the classi- 
cal ones with one electrón in the center of a pentagon. 
The quantum effects make the most probable coordinates 
slightly larger than the classical ones, but the difference 
vanishes for weak confinement. One can also see that p is 
more localized for weaker confinement, but only the one 
with ñu! = 0.01 meV looks like a Wigner molecule. In 
that case, all the points on the line p = 0.01 are closer 
to r* than any other r*. Also the snapshot animations 



of the Systems during the simulationsEJ show the case 
oí hu! = 0.01 meV to resemble what one expects for a 
Wigner molecule. As we work with a finite system, we 
do not have a real phase transition between Tilo = 0.1 and 
0.01 meV. In a two-dimensional electrón gas, the phase 
transition to a Wigner crystaUias been suggested to hap- 
pen at the density « 37±5a (in units of effective Bohr 
radius a|j, which is around 9.79 nm for our system param- 
eters). The radius of a circle that endoses one electrón 
on the average is thus around 360 nm at the transition 
point. This is in a good agreement with the inter-electron 
distances shown in Fig. ^, since in Fig. ^(c) the distance 
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FIG. 4: Conditional probability density p(r) for the right-most electrón. The contours are uniformly from 0.01 to 0.91. We 
mark with a plus the most probable electrón positions, and with a circle the classical positions. The confinement strength hu; 
is: (a) 10 meV, (b) 1 meV, (c) 0.1 meV, and (d) 0.01 meV. 



of edge-electrons to the center one is smaller than the 
suggested critical 360 nm, and in Fig. ^(d) the distance 
is roughly four times the critical one. The approxima- 
tive relation between and huj presented in Ref. ^ gives 
hoj « 0.034 meV for — 37, which is nicely between the 
cases of Figs. ||(c) and (d). 

It is also interesting to see the similarity between the 
efFective single-particle potentials V shown in Fig. ||(a) 
and (b), and the corresponding p in Fig. ||(c) and (d). 
This similarity could be uscd to describe the system, 
as we did above, as a system of six independent elec- 
trons, each with its own effective potential V. The as- 
sumption of parabolic V could easily be replaced by, e.g., 
V = ^{Cjl^x^+üjy^^yf) keeping the problem still solvable. 
There are, however, contributions that are not taken into 
account in this simple model, such as the exchange energy 
which is important for the spin-polarization of the sys- 
tem. Another aspect of the similarity between p and V is 
that it draws a connection between the wave function and 
the potential in a similar fashion that is often assumed 
in a semi-classical approximation. This could be used to 
motívate the studies of classical charged partióles in vari- 
ous confinements and also a semiclassical approach in the 



limit where the electrons are cióse to forming a Wigner 
molecule. 

The transition to a Wigner-molecule-like regime can 
also be seen in the radial pair distribution function, 
shown in Fig. ^. The function is clearly more peaked for 
weaker confinement. The peak at r = 1 consists of two 
types of electrón pairs, namely ones with both electrons 
on the edge, and ones with an electrón in the center and 
other on the edge. The electrón-electrón distance in the 
pair of the first type is in a classical solution 18% longer 
than in the second, the number of different pairs being 
the same in the classical solution. This double-nature 
cannot be seen in the first peak of g{r). This is not sur- 
prising, afHparticle exchanges happen even with huj = 
0.01 meVJl3 To study exchange, we have foUowed the 
most probable electrón positions while forcing the elec- 
trón originally in the center to move to the edge. The 
symmetry is broken by making small random displace- 
ments of the electrons. The particle exchange involves 
three or more electrons, as shown in Fig. One can 
see a collective rotation of the edge electrons before one 
is moved halfway to the edge. This kind of exchange is 
very easy for the small QDs as the one in question. One 
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FIG. 5: Potentials felt by one of the electrons (contours uni- 
formly from half to three) when the five other electrons are 
on their classical positions (marked with a circle). The con- 
finement strength huí is: (a) 0.1 meV, and (b) 0.01 meV. 



could argüe that larger electrón numbers would make the 
multi-particle exchanges to a more unlikely process. 

Compared with the experimental realizations of GaAs 
QDs, the electrón density where the Wigner molecule 
is found is extremely small. We feel that impurities 
would move the transition to larger densities as in the 
2D electrón ga&_where the transition is found to mov/ 
to Ts ~ 7.5 a^.ll3 For this, the approximative relatio: 
gives hu 0.37 meV, which is already closer to the typ- 
ical confinement strengths in experiments. It would also 
be very interesting to study the effect of impurities on 
the spin-polarization transition. 



IV. SUMMARY AND CONCLUSIONS 

We have first shown that the VMC method results 
energies in good agreement with the most accurate re- 
sults available. In VMC, the construction of the wave 
function clearly plays a central role. We have shown 
that the efficiency of SGA allows us to carefuUy optimize 
also the single-particle part of the wave function, result- 




FIG. 6: Radial pair distribution function for ñuj =1.0, 0.1, 
0.01 meV. The curves with smaller huj are more peaked. r is 
scaled to set the first peak to one. 




FIG. 7: Partióle exchange for huj = 0.01 meV. The electrón 
from the center (marked with light-gray circle) is moved to 
the right along the solid line. The rest of the electrons are 
always at their most probable positions. For the electrons 
ending at the starting position, only these positions and the 
path foUowed is shown. The circles showing the distance from 
the center have a radius from 500 to 1500 nm. 



ing VMC energies more accurate than the previous ones 
where the,aingle-particle states are taken from mean-field 
approach.EI On the other hand, our results show that in 
many cases the non-interacting single-particle states are 
optimal or very cióse to the optimal ones. This is proba- 
bly related to the high symmetry of the parabolic QD and 
the separation of the center-of-mass and relative motion. 
We feel that the efñciency of the SGA method would be 
even more useful in low-symmetry dots. 

Unlike in the previous accurate study of the six- 
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electrón QDp we have been able to reach low enough 
densities to find a spin polarization of electrons and, in 
an even lower density, a smooth transition to a Wigner- 
molecule-like state. The transition happeps roughly at 
the same density as in the 2D electron-gas.llá One should 
note that the 2D elcctron-gas does not spin-polarize be- 
fore the transition to a Wigner-crystaL_,but the spin- 
polarized state is very cióse in energy.Ej One possible 
explanation for the difference could be that the multi- 
particle exchange we find in the six-electron dot favors 
spin-polarization. In the 2D electron-gas, such process 
is less favorable. To qualitatively study the transition to 
a Wigner-molecule-like state, we have introduced a mea- 
sure function p, which we show to be very useful for the 
study of the electrón localization. 



Overall, we have found VMC to be a perfect tool for 
studying the properties of QDs in a wide range of system 
parameters, resulting energies in good agreement with 
the most accurate results available and enabling us to 
study the delicate transition of a QD to the classical 
regime. 
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